home *** CD-ROM | disk | FTP | other *** search
/ Software Vault: The Gold Collection / Software Vault - The Gold Collection (American Databankers) (1993).ISO / cdr48 / bpl70n12.zip / ARISOURC.ZIP / FPLOG.ASM < prev    next >
Assembly Source File  |  1993-03-07  |  7KB  |  172 lines
  1.  
  2. ; *******************************************************
  3. ; *                                                     *
  4. ; *     Turbo Pascal Runtime Library Version 7.0        *
  5. ; *     Real Logarithm                                  *
  6. ; *                                                     *
  7. ; *     Copyright (C) 1989-1993 Norbert Juffa           *
  8. ; *                                                     *
  9. ; *******************************************************
  10.  
  11.              TITLE   FPLOG
  12.  
  13.  
  14. CODE         SEGMENT BYTE PUBLIC
  15.  
  16.              ASSUME  CS:CODE
  17.  
  18. ; Externals
  19.  
  20.              EXTRN   RealAdd:NEAR,CmpMantissa:NEAR,RealFloat:NEAR,RealSub:NEAR
  21.              EXTRN   RealDivRev:NEAR,RealMulNoChk:NEAR,RealPoly:NEAR
  22.              EXTRN   HaltError:NEAR,ROverflow:NEAR,realmulfnochk:near
  23.              EXTRN   ShortMulRev:NEAR
  24. ; Publics
  25.  
  26.              PUBLIC  RLn
  27.  
  28.              IFDEF   EXTENSIONS
  29.              PUBLIC  RLog2,RLog10
  30.              ENDIF
  31.  
  32. ;-------------------------------------------------------------------------------
  33. ; RLn computes the natural logarithm of its argument. It uses a polynomial
  34. ; approximation to compute the natural logarithm of the reduced argument z. The
  35. ; reduced argument satisfies the inequality |z| <= (sqrt(2)-1)^2. RLog10 and
  36. ; RLog2 are additional routines that compute the logarithms base two and ten,
  37. ; respectively. Both first execute RLn to compute the natural logarithm and
  38. ; then proceed to multiply the result with the appropriate constants to get
  39. ; Log10 and Log2. The following polynomial approximation is used to compute
  40. ; the natural logarithm:
  41. ;
  42. ; rz := ((((0.09790802001953*z^2 + 0.1108818338371)*z^2 + 0.1428605246897)*z^2
  43. ;           0.1999999783036)*z^2 + 0.3333333333786)*z^2 * z + z
  44. ;
  45. ; This approximation has a theoretical maximum relative error of 3.20e-14.
  46. ; Maximum observed error when evaluated in REAL arithmetic is 9.31e-13.
  47. ;
  48. ; If the argument is negative or zero, runtime error 207 is invoked through the
  49. ; error handler.
  50. ;
  51. ; INPUT:     DX:BX:AX  argument
  52. ;
  53. ; OUTPUT:    DX:BX:AX  ln, log10, log2 of argument depending on routine called
  54. ;
  55. ; DESTROYS:  AX,BX,CX,DX,SI,DI,Flags
  56. ;-------------------------------------------------------------------------------
  57.  
  58.              IFDEF   EXTENSIONS
  59.  
  60. RLog10       PROC    FAR
  61.              MOV     DI,OFFSET $log_ten; push address of log10 tail-routine
  62.              JMPS    $start_log        ; compute common logarithm
  63. RLog10       ENDP
  64.  
  65.              ALIGN   4
  66.  
  67. RLog2        PROC    FAR
  68.              MOV     DI,OFFSET $log_two; push address of log2 tail-routine
  69.              JMPS    $start_log        ; compute logarithm dualis
  70. RLog2        ENDP
  71.  
  72.              ENDIF
  73.  
  74.              ALIGN   4
  75.  
  76. RLn          PROC    FAR
  77.              MOV     DI,OFFSET $log_end; push address of ln tail-routine
  78. $start_log:  OR      DH, DH            ; x negative ?
  79.              JS      $range_err        ; yes, error
  80.              OR      AL, AL            ; x zero ?
  81.              JZ      $range_err        ; yes, error
  82.              PUSH    DI                ; save log routine tail address
  83.              MOV     CX, 0FA81h        ; CL = exponent of constant a = 1,
  84.              MOV     SI, 0F333h        ;  DI:SI:CH = mantissa
  85.              MOV     DI, 03504h        ;   of 0.5*sqrt(2)
  86.              CALL    CmpMantissa       ; compare mantissas of x and 0.5*sqrt(2)
  87.              JNC     $gt_root2         ; if mantissa x > mantissa 0.5*sqrt(2)
  88.              DEC     CX                ; exponent of constant a = 0.5
  89.              DEC     AX                ; exponent = exponent - 1
  90. $gt_root2:   PUSH    AX                ; save exponent of x
  91.              MOV     AL, 80h           ; x = mantissa of x
  92.              XOR     CH, CH            ; clear LSB of constant a
  93.              PUSH    CX                ; save exponent of constant a
  94.              XOR     SI, SI            ; real constant
  95.              MOV     DI, SI            ;  a = 1 or a = 0.5
  96.              CALL    RealSub           ; x-a
  97.              POP     CX                ; get exponent of constant a
  98.              PUSH    DX                ; save
  99.              PUSH    BX                ;  x-a
  100.              PUSH    AX                ;   on stack
  101.              INC     CX                ; create
  102.              XOR     SI, SI            ;  constant
  103.              MOV     DI, SI            ;   2a
  104.              CALL    RealAdd           ; compute (x-a) + 2a = x+a
  105.              POP     CX                ; get
  106.              POP     SI                ;  back
  107.              POP     DI                ;   x-a
  108.              CALL    RealDivRev        ; compute (x-a)/(x+a)
  109.              MOV     CX, 5             ; polynomial has five coefficients
  110.              MOV     DI,OFFSET LN_COEFF; pointer to first coefficient
  111.              XOR     SI, SI            ; polynomial of type P(x^2)*x+x
  112.              CALL    RealPoly          ; z+z*p(z^2), max. rel. err. 2.6e-12
  113.              ADD     AL, 0FFh          ; compute rz := 2 * (z + z * p(^2))
  114.              ADC     AL, 1             ;  except when result is zero
  115.              POP     CX                ; get exponent
  116.              PUSH    DX                ; save
  117.              PUSH    BX                ;  rz on
  118.              PUSH    AX                ;   stack
  119.              XCHG    AX, CX            ; AL = exponent
  120.              SUB     AL, 80h           ; compute n = exponent - $80
  121.              CBW                       ; convert n to word
  122.              CWD                       ; convert n to longint
  123.              CALL    RealFloat         ; compute float (n)
  124.              MOV     CX, 0D280h        ; load
  125.              MOV     SI, 017F7h        ;  real constant
  126.              MOV     DI, 03172h        ;   ln(2)
  127.              CALL    ShortMulRev       ; compute n*ln(2),max. rel. err. 1.12e-12
  128.              POP     CX                ; get
  129.              POP     SI                ;  rz from
  130.              POP     DI                ;   stack
  131.              JMP     RealAdd           ; compute rz + n * ln(2)
  132.  
  133.              IFDEF   NOOVERFLOW
  134.  
  135. $range_err:  MOV     CH, -1            ; result negativ
  136.              JMP     ROverflow         ; largest REAL number
  137.  
  138.              ELSE
  139.  
  140. $range_err:  MOV     AX, 0CFh          ; load error code 207
  141.              JMP     HaltError         ; execute error handler
  142.  
  143.              ENDIF
  144.  
  145.              IFDEF   EXTENSIONS
  146. $log_ten:    MOV     CX, 0377Fh        ; load
  147.              MOV     SI, 0D8A9h        ;  constant
  148.              MOV     DI, 05E5Bh        ;   1/ln(10)
  149.              JMPS    $mult_const       ; compute common log from natural log
  150. $log_two:    MOV     CX, 05C81h        ; load
  151.              MOV     SI, 03B29h        ;  constant
  152.              MOV     DI, 038AAh        ;   1/ln(2)
  153. $mult_const: CALL    RealMulNoChk      ; compute log dualis from natural log
  154.              ENDIF
  155.  
  156.              ALIGN   4
  157.  
  158. $log_end:    RET                       ; done
  159.  
  160. LN_COEFF     DB      07Dh,               084h,048h  ;  9.790802001953e-2
  161.              DB      07Dh,068h,0D0h,003h,016h,063h  ;  1.108818338371e-1
  162.              DB      07Eh,0BAh,085h,007h,04Ah,012h  ;  1.428605246897e-1
  163.              DB      07Eh,00Fh,058h,0CBh,0CCh,04Ch  ;  1.999999783036e-1
  164.              DB      07Fh,00Eh,0ABh,0AAh,0AAh,02Ah  ;  3.333333333786e-1
  165. RLn          ENDP
  166.  
  167.              ALIGN   4
  168.  
  169. CODE         ENDS
  170.  
  171.              END
  172.